import numpy as np
import pandas as pd
import math
from scipy import stats
from scipy.stats import norm, skew #for some statistics
import matplotlib.pyplot as plt
%matplotlib inline
import seaborn as sns
sns.set()
from sklearn.model_selection import train_test_split
from sklearn.linear_model import LinearRegression
from sklearn.metrics import mean_squared_error, r2_score
import lightgbm as lgb
### Functions =================================================================
def xdescribe(df, dsp=True, ret=False): # pandas.describeメソッドの拡張
print(f'## Shape:{df.shape}')
dtype = pd.DataFrame({'dtype' : df.dtypes })
n_na = pd.DataFrame({'NAs' : df.isnull().sum()})
n_unique = pd.DataFrame({'n_unique': [len(df[c].unique()) for c in df.columns]},
index=df.columns)
df_descr = df.describe(include='all').T
if 'unique' in df_descr.columns:
df_descr = df_descr.drop(['unique'], axis=1)
df_descr = pd.concat([dtype, n_na, n_unique, df_descr], axis=1)
df_descr['count'] = df_descr['count'].astype(int)
if dsp:
# backup & restore sttting
mr = pd.get_option('display.max_rows') # backup sttting
pd.set_option('display.max_rows', None) # unlimitted
display(df_descr)
pd.set_option('display.max_rows', mr) # restore sttting
prt = True
for idx in df_descr.index:
if df_descr['dtype'].astype(str)[idx] == 'category':
if prt:
print(f'### Featur: Categories')
prt = False
print(f' {idx}: \t{df[idx].cat.categories}')
if ret:
return df_descr
else:
return
### Functions =================================================================
# 目的変数の分布を確認
def target_dist(df, target_colname):
# Get the fitted parameters used by the function
(mu, sigma) = norm.fit(df[target_colname])
skewness = skew(df[target_colname].dropna())
print( f'\n{target_colname}; 平均={mu:.2f}, 分散={sigma:.2f}, 歪度={skewness:.2f}')
# plot the distribution
fig = plt.figure(figsize=(12,4))
ax = fig.add_subplot(1, 2, 1)
ax.set_title(target_colname+' distribution')
sns.histplot(df[target_colname], stat='density', bins=20,
kde=True, line_kws={'linewidth':2, 'color':'black'})
t_min, t_max = df[target_colname].min(), df[target_colname].max()
sp = (t_max-t_min)/100
x = np.arange(t_min-sp*10, t_max+sp*10, sp)
y = norm.pdf(x, mu, sigma)
sns.lineplot(x=x, y=y, color='black')
# Q-Q plot
ax = fig.add_subplot(1, 2, 2)
stats.probplot(df[target_colname], plot=plt)
plt.show()
### Functions =================================================================
# 目的変数(target)が連続、説明変数(Feature)がカテゴリ変数変数の場合
# カテゴリ変数の分布を確認
def plt_box_dist(df, target, feature):
n_unique = len(df[feature].unique())
font_size = 12 - n_unique / 5
fig_height = 1 + n_unique / 4
fig = plt.figure(figsize=(8,fig_height))
# boxplot
ax0 = plt.subplot2grid((1, 6), (0, 0), colspan=4)
sns.boxplot(y=df[feature], x=df[target], orient='h', color='#628CD2', width=0.7,
flierprops={'markersize':5,'markerfacecolor':'gray',
'markeredgecolor':'None',
'alpha':min(1, 1/np.log10(df.shape[0]/10+9))})
ax0.set_yticklabels(ax0.get_yticklabels(), size=font_size)
if pd.api.types.is_numeric_dtype(df[feature]):
ax0.invert_yaxis()
# distribution
ax1 = plt.subplot2grid((1, 6), (0, 4), colspan=2)
sns.countplot(y=df[feature], color='#628CD2')
ax1.set_yticklabels([])
ax1.set_ylabel('')
if pd.api.types.is_numeric_dtype(df[feature]): # 数値の場合には下から
ax1.invert_yaxis()
plt.show()
### Functions =================================================================
# 目的変数(target)が連続、説明変数(Feature)が連続(数値)変数の場合
# 目的関数との相関/分布をプロット
def plt_reg_dist(df, target, feature):
n_unique = len(df[feature].unique())
font_size = 10
fig_height = 3
fig = plt.figure(figsize=(8,fig_height))
# regplot
ax0 = plt.subplot2grid((1, 6), (0, 0), colspan=4)
sns.regplot(data=df, x=target, y=feature, color='#3050A0',
scatter_kws={'alpha':1/np.log10(df.shape[0]+9), 's':12},
line_kws={'lw':1, 'color':'black'})
# distribution
ax1 = plt.subplot2grid((1, 6), (0, 4), colspan=2)
sns.histplot(data=df, y=feature, stat='count', bins=30, color='#628CD2',
shrink=1.0, kde=True)
ax1.set_yticklabels([])
ax1.set_ylabel('')
plt.show()
### Function =================================================================
# 予測誤差(true-pred)プロット
def true_pred(df, true_colname, pred_colname):
fig, ax = plt.subplots(figsize=(5,5))
xymin = min(df[true_colname].min(), df[pred_colname].min())
xymax = max(df[true_colname].max(), df[pred_colname].max())
pitch = (xymax-xymin)/100
sns.scatterplot(data=df, x=true_colname, y=pred_colname,
color='mediumblue', alpha=1/np.log10(df.shape[0]+9))
ax.set_xlabel(true_colname)
ax.set_ylabel(pred_colname)
ax.set_xlim([xymin-pitch*10, xymax+pitch*10])
ax.set_ylim([xymin-pitch*10, xymax+pitch*10])
x = np.arange(xymin-pitch* 10, xymax+pitch*10, pitch)
sns.lineplot(x=x, y=x, color='red', linewidth=1, linestyle='--')
x = np.arange(xymin-pitch* 5, xymax+pitch* 5, pitch)
reg = LinearRegression().fit(df[[true_colname]], df[pred_colname])
sns.lineplot(x=x, y=reg.coef_[0]*x+reg.intercept_, color='red', linewidth=2)
RMSE = np.sqrt(mean_squared_error(df[true_colname], df[pred_colname]))
R2 = r2_score(df[true_colname], df[pred_colname])
plt.text(x=0.1, y=0.85, transform=ax.transAxes, fontsize=12,
s=f'RMSE = {RMSE:.3f}\nR2 = {R2:.3f}')
### ===========================================================================
### ===========================================================================
### Main
df = sns.load_dataset('tips') # seaborn tipsデータセット
df_xdesc = xdescribe(df, ret=True) # 各変数の属性を表示
target_colname = 'tip' # 目的変数を指定
print(f'\n## 目的変数の分布を確認')
target_dist(df, target_colname)
print(f'\n## 説明変数の分布を確認')
feature_colnames = [c for c in df.columns if c !=target_colname]
for i, colname in enumerate(feature_colnames):
print(f'\n## {i+1}/{len(feature_colnames)}, {colname}')
if str(df_xdesc.loc[colname, 'dtype']) in ['object', 'category']:
plt_box_dist(df, target_colname, colname)
else:
plt_reg_dist(df, target_colname, colname)
if df_xdesc.loc[colname, 'n_unique'] <= 30:
plt_box_dist(df, target_colname, colname)
# regression
# データフレームを学習・評価、説明変数(特徴量)と目的変数に分割
df_train, df_test = train_test_split(df, train_size=0.75, random_state=123)
X_train = df_train.drop(target_colname, axis = 1)
y_train = df_train[target_colname]
X_test = df_test.drop(target_colname, axis = 1)
y_test = df_test[target_colname]
lgb_train = lgb.Dataset(X_train, y_train) # 学習用
lgb_test = lgb.Dataset(X_test, y_test) # 評価用
# パラメータの設定
params = {'objective':'regression', 'metric':['rmse', 'r2'], 'verbosity':-1}
# LightGBMで学習
lgb_model = lgb.train(params, lgb_train)
# 評価データセットで予測
y_pred = lgb_model.predict(X_test)
df_pred = pd.DataFrame({'true': y_test, 'predicted': y_pred})
print('## LightGBM regression結果')
# 予測誤差(true-pred)プロット
true_pred(df_pred, 'true', 'predicted')
